A multi-region nonlinear age-size structured 
fish population model 



Blaise Faugeras * and Olivier Maury 

IRD, CRHMT, Av. Jean Monnet, BP 171, 34200 Site, France 

Abstract 

The goal of this paper is to present a generic multi-region nonlinear age-size struc- 
tured fish population model, and to assess its mathematical well-posedness. An 
initial-boundary value problem is formulated. Existence and uniqueness of a positive 
weak solution is proved. Eventually, a comparison result is derived : the population 
of all regions decreases as the mortality rate increases in at least one region. 
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1 Introduction 



Fish population dynamics models are essential to provide assessment of the 
fish abundance and fishing pressure. Their use forms the basis of scientific 
advice for fisheries managements. Discrete age structured models are most of 
the time used for fisheries stock assessments [1] . Indeed, ecologists, mathemati- 
cians and population biologists have observed that the age structure provides 
more realistic results at reasonable computational expense for a wide variety 
of biological populations (see [2], [3], [4], [5]). 

In this paper we study a model which was first designed to represent Atlantic 
bigeye tuna populations [6] but which is also generic enough to be potentially 
usefull for various fish species. Indeed, most fish populations share specific 
characteristics which have to be taken into account in order to model their 
dynamics in a realistic manner. 

A first point concerning tuna fisheries is that they are highly heterogeneous 



Email address: Blaise.Faugeras@ifremer.fr (Blaise Faugeras). 



Preprint submitted to Elsevier Science 



1 July 2010 



in space and time. This has an important impact on their functioning. Impor- 
tant migrations of fish occur at various scales and fish movements have to be 
explicitly represented. Moreover, growth potentially varies with space that is 
to say with the region of the ocean under consideration. Hence, fishes of the 
same age can exhibit very different sizes depending of their various history. 
Consequently a spatialized approach taking explicitly into account the poten- 
tial variability of growth in space has to be used. 

A second point is that, because of non- uniform mortality over sizes, bias on 
both growth and mortality estimates may result from simply adding a gaus- 
sian size distribution to an age structured model as it is generally done. It is 
reasonable to think that the use of both age and size as structure variables 
should enable to overcome this difficulty 

These are some of the principal problems of current stock assessment models. 
That is why it is necessary to carry on the modelling effort by proposing and 
testing more complex models. This paper follows this direction and its purpose 
is twofold. 

First we describe a synthetic and generic model of population dynamics in 
which both age and size are taken as structure variables and in which fish 
movements among spatial regions are explicitly represented. The model is a 
system of coupled partial differential equations. Nonlocal nonlinearities appear 
in the boundary conditions modelling recruitment that is to say the birth law 
or density dependent fish reproduction. The relative complexity of the model 
enables a direct and simultaneous comparison with all the data available for 
tuna fisheries such as catches, fishing efforts, size frequencies, tagging data, 
and otoliths increments. This paper does not aim at getting into all the details 
of the parameterizations used to represent a particular tuna population and 
we refer to [6] for these points. 

Our second and most important goal is to assess the mathematical well- 
posedness of the model. The paper is organized as follows. The equations of 
the model are presented in Section 2. Sections 3, 4 and 5 deal with the math- 
ematical analysis of the model. In Section 3 we formulate an initial-boundary 
value problem, introduce a variational formulation and state our main math- 
ematical results. Existence of a unique weak solution is shown in Section 4. 
As often with nonlinear problems the proof uses a fixed point argument. The 
methodology follows the one proposed in [7] for a scalar equation. It has to be 
adapted in order to be able to deal with our nonlinear system. We also show 
positivity of the solution and give a comparison result in Section 5. Namely 
we prove that if the fish mortality rate increases in at least one geographic 
region then the population globally decreases in all regions. 
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2 The model 



The dynamics of the population of fish is described through density functions 
Pi(t,a,l) where time t G (0,T), age a G (0, A) and length I G (0, L) are 
continuous variables and the subscript i G [1 : N] refers to the geographic 
zone or region under consideration. The number of fish of age between oi and 
a2, of length between l\ and l 2 at time t in region i is given by the integral 



Let us set O = (0, T) x (0, A) and Q = O x (0, L). The time evolution of the 
population given by Eq. 1 includes the following processes. 

In region i, as time goes on and fishes grow older, their length increases with a 
growth rate 7$. In a fish population individuals of the same age can often differ 
markedly in size [8] . This variability in growth can result from many different 
mechanisms, including genetic or behavorial traits that confer different per- 
formances to individuals, and factors such as environmental heterogeneity and 
variability [9]. In fishery science, this variability is usually taken into account 
in age-structured models using a length-at-age relation perturbed by a Gaus- 
sian noise with a length dependent standard deviation (see for example [10]). 
The model discussed here is length-structured and uses a diffusion term in the 
length variable with dispersion rate di to account for individuals having the 
same age but different lengths. The advection-diffusion term in length can be 
seen as the limit of a random walk model in which each individual grows with 
an average velocity, but has at each time step a small binomial probability to 
grow faster or slower than this average (see the book by Okubo [11] for more 
details). 

The model also describes mortality and migration of individuals. The mortality 
rate is split into natural mortality Hi and fishing mortality /«. Let also 
be the migration rate of individuals going from region i to region j {m^j = 
if regions i and j are not adjacent). 

The density functions pi for i G [1 : N] follow the balance law: 




d t pi(t,a,l) +d a pi(t,a,l) 



di(di(t, a, l)dipi(t, a, I)) - di(^i(t, a, l)pi(t, a, I)) 



N 



+ ^2m j ^ i (t,a,l)p j (t,a,l) 



< 
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These equations have to be completed with initial and boundary conditions. 
Homogeneous Neumann boundary conditions at I — and I = L express the 
fact that the length of individuals can not reach negative values or values 
larger than L. 

d lPi (t,a,0) = d lPi (t,a,L) = 0, (t,a)eO. (2) 

The initial age and size distribution is prescribed, 

Pt (0, a, I) = p°(a, I), (a, I) e (0, A) x (0, L). (3) 

We also need a boundary condition for a = that is to say a recruitment law. 
It is written as: 

Pi (t, 0, 1) = fat, I, Pi(t)), (t, I) e (0, T) x (0, L), (4) 

The length of recruited fish is assumed to lie between and a small constant 
length Lb. Moreover we denote by L m the minimal length of fishes which have 
reached maturity. and L m satisfy < Lb < L m < L. The stock spawning 
biomass is calculated as 

p i(t) = [ f Wi(t,a,l)pi(t,a,l)dlda, (5) 

JO JL m 

where iUj is a weighting function. Finally we use a Beverton and Holt [9] 
stock-recruitment relation in each region and obtain, 

A(t,/,F)=l [ o, L6] (/)^(t)^p, (6) 

wherel[ 0i L 6 ] is the usual characteristic function, Q { > is a constant parameter 
and ipi(t) is a given function of time used to parameterize fluctuations of the 
recruitment not taken into account in the Beverton and Holt relation. 



3 Main assumptions and preliminary results 

In this section we set the mathematical frame in which the analysis is con- 
ducted. We formulate the main assumptions which are made on the data of 
the model, give the definition of a weak solution to the initial-boundary value 
problem and state our results in Theorems 3.1 and 3.2. 

3. 1 Functional spaces 

Let us introduce the functional spaces which we use in the remainder of this 
work. 
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The vectorial notation p = (pi, ...,pn) t is used. The usual scalar product of 
two vectors p, q G M. N is denoted by p.q and the norm of p by |p|. 
H and H 1 are the separable Hilbert spaces defined by H = (L 2 (0, L)) N and 
H 1 = (-£^(0, L)) N . H is equipped with the scalar product 



(P,q)n= [ L p(l).q(l)dl. 
Jo 



We denote by ||.||h the induced norm on H. 
H 1 is equipped with the scalar product 

(p,q)w = f L p(l).q(l)dl+ ( L d^.d^dl. 
Jo Jo 

We denote by H-Hh 1 ; the induced norm on H 1 . 

By < ., . > we denote the duality between H 1 and its dual (H 1 )'. 

L 2 (C, H) (resp. L 2 (0, H 1 )) denotes the Hilbert space of measurable functions 

of O with values in H (resp. H 1 ) such that 

||pIU 2 (o,h) = (/ ||p(*,a, -)\\udtda) 1/2 < oo (resp. ||p||z,2(o,Hi) < oo). 

We also make use of the notation V = L 2 {0,'R l ) and the dual space V = 
L 2 (0, (H 1 )'). By << ., . >> we denote the duality between V and its dual 
V. 

The partial derivatives d t and d a denote differentiation in V(0, (H 1 )') and D 
stands for d t + d a . 

We will have to use he following trace result. 

Lemma 3.1 Let p, q e V such that Dp, Dq e V . It holds that: 
For all t G (0,T) and all a G (0,A), p has a trace at t — to belonging to 
(L 2 ((0,A) x (0,L))) N and at a = a belonging to (L 2 ((0,T) x (0,L))) N . The 
trace applications are continuous in the strong and weak topology. Moreover 
the following integration by parts formula holds, 



( [< £>p,q > + < Dq, p >]dtda = [ [ [p.q(T, a, I) - p.q(0, a, l)]dadl 
Jo Jo Jo 

+ [ T f L [p.q(t,A,l) -p.q(t,0,l)}dtdl 
Jo Jo 



Proof: This result is the extension to dimension iV of Lemma in [12]. Also 
see [13]. □ 



We will also have to consider the space L°° = (L°°(Q)) N . L°°(Q) is a Banach 
space equipped with the norm ||pi||oo — inf{M; \pi(t, a, l)\ < M a.e. in Q}. 
Similarly L°° is a Banach space equipped with the norm ||p||oo = max ||pj||oo- 

i£[l:N] 
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3.2 Assumptions on the data and preliminary transformation of the system 



The movements rates m^j are assumed to satisfy 

• m i^-j{t, a, I) > a.e in Q, G 

We define the matrix of movements M by 



M i:j = { 



rrij^i if i ^ j, 

N 

- m i~>k if i = j. 



N N 

Hence the term, rrij^iPj — m^)^], in Eq. 1 can be written in matrix 
form as (Mp) ; . 

Concerning the diffusion coefficients dj, the growth rates 7, and the natural 
and fishing mortality rates fa and fi, we make the following assumptions for 
all i G [1 : N]: 

• di{t,a,l) >d°>0, a.e in Q, d, t G L°°{Q), 

• 7j(t,a, I) is differentiable with respect to /, and 7i,<9;7i G L°°(Q), 

• fa(t,a,l), fi(t,a,l) > 0, a.e in Q, fa, fi G L°°(Q). We also make use of the 
notation zt = fa + fi. 

In the formulation of the recruitment process (Eq. 4) and Wi satisfy: 

• i/>i(t) > a.e in (0,T) and ^ G L°°(0,T), 

• Wj(t, a, /) > a.e in Q and Wi G L°°(Q). 

The initial distributions p®(a,l) satisfies for all i G [1 : AT]: 

• p?(o, > a.e in Q, p° G L 2 ((0, A) x (0, L)). 

In order to prove our existence result it is convenient to perform a change of 
unknown function: p satisfies (l)-(4) if and only if p = e~ xt p is a solution to 
the same system where —(fa + fi)pi is replaced —(fa + fi + \)p% in Eq. 1 and 
fa in the expression of the boundary condition at a = (Eq. 4) is replaced by 

m^m)) = i [0 , Lb] m(t) - * {t \ uV (7) 

Vie-™ + Pi(t) 

Pi(t) — / / Wi(t,a : l)pi(t,a : l)dlda. (8) 

Jo JL m 

In the remaining part of this paper this change of unknown is implicitly done 
and we omit the pi notation. The constant A will be fixed to a convenient value 
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below. Moreover, the possible nullification of the term O^e xt + Pi(t), invites 
us to define, 

P(t) 

mi,m)) =v>i^(t) g, e _ M ;| p , (f)| . (9) 

This formulation will be used in the following. We will show that if initial dis- 
tributions, p® are nonnegative then Pi > a.e. in Q, thus the two formulations 
are equivalent. 



3. 3 Variational formulation and weak solutions 

Formally multiplying Eq. 1 by a function qi and integrating by parts on (0, L) 
yields to the definition of the following linear forms. For Pi,qi G H l (0, L) let 
us define, 

rL rL rL 

bi(Pi,Qi)= didipidiqidl + li{diPi)qidl + j + 8^ + X)piqidl, (10) 

JO JO JO 

Ci(p,ft) = - [ L (Mp) iqi dl, (11) 

JO 

e»(P, Qi) = h{Ph Qi) + Ci(p, qi) (12) 
Summing over i, we define for p, q G H 1 , the bilinear form e(p, q) by, 

N 

e(p,q) = 5Ze i (p,g i ) (13) 

i=i 

Lemma 3.2 For X > (^jqIItIIoo + II^Tlloo + -^llMHoo), the bilinear form 

e(., .) is continuous and coercive on H 1 x H 1 , i.e there exist constants C± > 
and C2 > such that 



|e(p,q)| < Ci||p||Hi||q||Hi, Vp.qGH 1 , (14) 
e(p,p) >C 2 ||p||^, VpGH 1 . (15) 
Proof: Using Cauchy-Schwarz inequality we obtain, 

\^2h(Pi,qi)\ < (||d||oo + | |t| |oo + INU + Hdrrl U + A)||p|| H i||q| 



and 

'o 



/■L 

^Q(p,gi)| = 1X1X1 / M ijPj1idl\ < | |]M| looA^I |p| | H i I |q| |ni • 
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which proves (14). 

Again using Cauchy-Schwarz inequality yields 

I / JidiPiPi\dl < ||7||oo||^Pj||L2(0,L)lbi||L2( ,L). 
J 

Young's inequality then gives for any a > 

rL q, | 

I J iidiPiPidl\ < -\\diPi\\ 2 L 2^ L) + —\\l\\lo\\Pi\\h(o,L)- 
Therefore we have that 



Y,J Q li d iPiPi dl > -f II%>IIh - 2^ll7llLllP 



2 

oollf IIH- 



Now since \ii and fi are positive and di is bounded below by d°, it follows that 
e(p,p) > (d° - |)Pp||^ + (A - (t^HL + \\da\L + ^|M|U))||P||^ 

It is possible to choose a = d° and A such that A = (A — (— -| I71 1^ + 1 \d^\ \ 2 OQ + 

ZiOj 

(f 

iVUMlloo)) > and C 2 = min( — ,X°). □ 



We can now give the definition of a weak solution to the initial-boundary value 
problem (l)-(4) and state the results which are shown in Section 4 and 5. A 
weak solution to the initial-boundary value problem (l)-(4) is a vector valued 
function p satisfying the following problem (P): 
Find 

p G V, such that Dp G V , (16) 

solution of 

/ < Dp, q > dtda + [ e(p, qjdtda = 0, Vq G V, (17) 
Jo Jo 

p(0,o,Z) = p°(a,Z) a.e in (0, A) x (0, L), (18) 
p(t,0,Z) = (3(t,l,P(t)) a.e in (0,T) x (0,L). (19) 

In Section 4 it is proved that: 

Theorem 3.1 There exists a unique solution p to problem (P). 

Notation: p(t, a, I) and q(t, a, I) being vector valued functions, p < q means 

that Pi < qi a.e. in Q for all i G [1 : N]. 

With this notation, it is proved in Section 5 that: 
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Theorem 3.2 The solution, p, to problem (P) is nonnegative a.e in Q. 
Moreover, let p 1 (resp. p 2 ) denote the solution to problem (P) associated 
with the vector of mortality rates z 1 (resp. z 2 ). If z 1 < z 2 then p 2 < p 1 . 



4 Existence and uniqueness 



The proof of existence and uniqueness consists in two main steps. First we 
show the result in the case of a constant recruitment (independent of the fish 
density). Second a fixed point argument enables to cope with the original 
nonlinear recruitment. 

Lemma 4.1 Let b be fixed in L 2 ((0,T) x (0,L)) N . There exists a unique p 
satisfying (16)- (18) of problem (P) in which the initial condition (19) is 
replaced by p(t, 0, 1) = b(t, I) a.e in (0, T) x (0, L). 

Proof: The proof is an adaptation of the results given for the scalar case in 
[7]. We sketch it for the sake of completeness. It consists in two steps. 

Step 1: We prove that given h G V there exists a unique p G V, Dp G V 
such that 

/ < Dp, q > dtda + / e(p, q)dtda = < h, q > dtda, Vq G V (20) 
Jo Jo Jo 

and p(0,o,Z) = p(t,0,l) = 0. 



Let A be the unbounded linear operator on (L 2 (Q)) N with domain D(A°) = 
{p E (L 2 (Q)) N , dtp + d a p E (L 2 (Q)) N , p(0, t, I) = p(t, 0, = 0}, defined by 
p G D(A°), A°p = d t p + d a p. Then — ^4° is the infinitesimal generator of a 
contraction semigroup, (S(r)p, t > 0), in (L 2 (Q)) N (see [14]) and 



(S(r)p)(t,a,l) 



p(t — r,a — t,1) if (t — r, a — r, I) G Q, 
otherwise. 



From this one can deduce that the unbounded linear operator A from V to V 
with domain D(A) = {p G V, Dp G V, p(0, a, I) = p(t, 0, 1) = 0}, defined by 
Ap = Dp is a maximal monotone operator. 

With the bilinear form e(., .) we can define a linear bounded and coercive op- 

erator E from V to V such that « Ep, q »= f e(p, q )dtda, Vp, q G V^. 

Jo 

Since is bounded and coercive and A is maximal monotone we conclude 
that for any h G V there exists a unique p G -D(A) solution to Ap + Ep = h 



9 



which is an abstract formulation of our problem because 
<< Ap, q >>= J < Dp, q > dtda, Vp G D(A), Vq G V. 

Step 2: Let us now introduce a sequence of functions 4> n G (C°°(Q)) Ar such 
that 

B (O, a, p°(a, in (L 2 ((0, A) x (0, L))) N , 

cf> n (t, 0, b(f , in (L 2 ((0, T) x (0, L))) N , 

From step 1, we conclude that there exists a unique q n in D(A) solution to 
Aq n + Eq n = -Acf> n - E<p n . Therefore p n = q n + n is a solution to (17) 
satisfying p n (0, a, I) = n (O, a, Z) and p n (t, 0, Z) = (j> n (t, 0, /). 
Now taking p n as a test function in (17), integrating by parts using Lemma 
3.1 and using the coercivity of e(., .) we obtain that 

C2||p n ||y < 2^ n ^' a '^tL 2 ((Q,A)x(0,L))) N + 2^ n ^ t,Q, ^t L2 ((0,T)><(0,L))) N - 

By the choice of 4> n this implies that p™ is a bounded sequence in V. There- 
fore we can extract a subsequence still denoted p™ such that p n — > p weakly 
in V and Dp n — > r weakly in V. Since the operator D is continuous on 
T>'(0, (H 1 )'), r = Dp. Moreover since E is continuous Ep n — > Ep. We con- 
clude that p satisfies (17). The continuity of the trace applications on t = 
and a = implies that p(0, a, I) = p°(a, I) and p(t, 0, 1) = b(t, I). □ 



Lemma 4.2 Let C 3 = ( max [AL 2 | |^| 1^1 |w| | 2 X) (—) 2 ]) 1/2 ; then the applica- 

ie[l:N] fj i 

tion 

(pi(t, a, I)) i — y (fii(t, I, Pi(t))) (cf Eqs 8 and 9) defines a bounded nonlinear oper- 
ator, lipschitz continuous from L 2 {0, H) to (L 2 ((0, T) x (0, L)) N with lipschitz 
constant C 3 . 

Proof: The application, pi(t,a,l) i-> P%(t) = / Wi(t,a,l)pi(t,a,l)dadl, 

Jo Jh m 

defines a bounded linear operator from L 2 (Q) to L 2 (0,T). This follows from, 

rA r-L r-A r-L 

| / / Wi(t, a,l)pi(t, a,l)dadl\ < / / \n)i(t, a, l)pi(t, a, l)\dadl, 
Jo JL m . Jo Jo 

and using Cauchy-Schwarz yields 

[ T \Pi(t)\ 2 dt < ||w||^L|| Pi ||| a(e) . 

J 

The application pi(t, a, I) >->■ j3i(t, I, Pi(t)) defines a bounded nonlinear operator 
from L 2 (Q) to L 2 ((0, T) x (0, L)). This follows from the fact that the applica- 
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P e XT 
tion Ui(t,P) = — — t — from [0,T] xRtoR satisfies k(t,P)| < — — 

and therefore we have 

J o j o (A(t,^^)) 2 ^<IIV'llL(^) 2 llw||^ 2 ll^lli 2 (Q)- 

Lipschitz continuity follows from the fact that (£, P) h- >■ m(t,P) is lipschitz 
continuous in P uniformly in £ G [0,T], 

AT 

k^P 1 ) - ^(t,P 2 )| < -HP 1 - P 2 |, VP\P 2 G K, V* G [0,T]. 
Hence, if to (resp. ) we associate P/ (resp. P?) it holds that 

T L 

f [ m,i,pi{t))-^{t^p^t))\ 2 dtdi 

Jo Jo 

^IMIL^) 2 / ip/w-p^wi 2 ^, 

&i Jo 

XT 

□ 



Lemma 4.3 TTiere exzsfo a unique p satisfying problem (P) 

Proof: Let p be given in V. With p we associate a vector (p(£)). Let us denote 
Pp = p the solution to (16)-(18) and satisfying (pi(t, 0, 1)) = (/3i(t, I, P(£))- 
From Lemma 4.1 and Lemma 4.2 we deduce that the nonlinear operator T 
maps V into itself. Moreover it follows from Lemma 3.1 that 

/ <Dp,p>dtda>~ [ ||p°(a,.)||^a-i F \\p(t,0, .)\\ 2 u dt 
Jo 2 Jo 2 Jo 

The coercivity of e(., .) leads to 

C 2 / \\p(t,a,.)\\ 2 ul dtda<\ [* \\p°(a, .)\\ 2 n da + \ f ||p(i, 0, .)\\ 2 u dt. 
Jo 2 Jo 2 Jo 

Lemma 4.2 then gives 

r 1 r A 1 

C 2 J G \\p(t,a,.)\\ 2 n idtda< -J ||p°(a, .)||h^o + t^HpII^^h), 

and JF is bounded from L 2 (0, H) to V. 

The solutions we are looking for are the fixed points of P. Let us show that 
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T is a strict contraction in L 2 (0,H). 

Let p 1 and p 2 be given in L 2 (0, H) and let p 1 = ^p 1 and p 2 = J-p 2 be 
the associated solutions. The difference p = J-p 1 — J-p 2 satisfies (16), (17), 
p(0, a, I) = and ( Pi (t, 0, 1)) = (fr(t, I, P*(t)) - I, Pf{t))). 
At the end of the proof of Lemma 3.2, since A is arbitrary, one can choose 

A = Ai + A 2 with Ai > ^olMliL + II^TllL + ^|| M IU and A 2 > arbitrary. 
Hence, 

e(p,p) > C 2 ||p||hi + A 2 ||p||h > A 2 ||p||h, Vp G H 1 . 
Now using Lemma 3.1 once again we obtain 



X 2 f \\p(t,a,.)\\ 2 n dtda<-Y, [ I m,l,PHt))-Pi(t,l,P?(t))] 2 dtdl 
Jo z ~[ Jo Jo 

and Lemma 4.2 gives 



A2||pIIl 2 (o,h) ^ 9C3IIP 1 - p 2 l 



L 2 (C,H)- 



We can choose A 2 = C 3 and since p = J-p 1 — Fp 2 this proves that J 7 is a 
strict contraction on L 2 (0,H). Thus it follows from Banach fixed point the- 
orem that T admits a unique fixed point p which is the desired solution. □ 



5 Positivity and comparison result 



In this section we first show in Lemma 5.1 that the fish density population 
solution to our model is positive. Then a comparaison result is given in Lemma 
5.2. 

Lemma 5.1 The solution p to problem (P) is nonnegative a.e. in Q. 

Proof: As in the proof of Lemma 4.3, let p be given in V and let J-p = p 
denote the solution to (16)-(18) and satisfying (pi(t,0,t)) = (Pi(t,l, Pj(t))). 
Let us also assume that p > 0. 

The negative parts of pi(0, a, I) and p { (t, 0, 1) satisfy (pi(0, a, = (p°(a, = 
Oand(p i (t,0,Z))- = (A(t,i,A(*)))- = 0. 

One can then show using Lemma 3.1 (see [7]) that / < Dp, p~ > dtda < 0. 

Jo 

The bilinear form e can be decomposed as e(p,p~) = e(p + ,p~) — e(p~,p~), 
with 

N 

e(p + ,P~) = J2 b i(Pt'Pi) +Ci(P + ,P7)- 

i=i 

It holds that h{pf ,p~) = since one can check that bi(pi,p~) = — bi(p~,p~). 
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Moreover, c i (p + ,p i ) = — / M^p\p~dl < 0, since M^p+pl > for i ^ j 

and Mupfp^ = 0. 

We conclude that e(p + ,p~) < 0. 

Taking q = p in Eq. 17 yields, 

J <Dp,p~>dtda + J e(p + , p~)dtda — J e(p~, p~)dtda = 

so that we obtain J* e(p _ , p~)dtda < 0. The coercivity of e gives, 
C^Hp Hv < 0, that is to say p is nonnegative. 

If we define a sequence with p 1 = p and p n+l = J-p n , then from the previous 
lines we deduce that p n is nonnegative for all n > 1. By Banach fixed point 
theorem this sequence converges to the solution p which is therefore nonneg- 
ative. □ 



Lemma 5.2 Let p 1 (resp. p 2 ) denote the solution to problem (P) associated 
with the vector of mortality rates z 1 (resp. z 2 ). If z 1 < z 2 then p 1 > p 2 . 

Proof: 

Step 1: Let p 1 and p 2 be given in V and satisfying < p 1 < p 2 . Let p 1 = J-p 1 
and p 2 = J-p 2 be the associated solutions defined as in the proof of Lemma 
5.1. Let us show that p 1 < p 2 . 

It is clear that Pl(t) < Pf{t) a.e in (0, T), then since Ui(t, P) is an increasing 
function of P it holds that p t (t, I, P±(t)) < I, P 2 {t)) a.e in (0, T) x (0, L). 
The difference p = p 2 — p 1 satisfies (16), (17) and 

P (o,u) = o, 

p(t, 0, a) = (frit, I, P 2 (t))) ~ (Pi(t, I, Pl(t))) > 0. 

This is the same situation as in the first part of proof of Lemma 5.1 and we 
conclude that p is nonnegative that is to say p 1 < p 2 . 

Step 2: Let p > be given in V. To the vectors of mortality rates z 1 and z 2 
(0 < z 1 < z 2 ) we associate the bilinear forms e 1 and e 2 (see (10)-(12), note 
that c] (.,.) = c 2 (.,.)) as well as the nonlinear operators J 71 and J 72 defined as 
in the proof of lemma 5.1. They define the solutions p 1 = F l p and p 2 = J^p. 
Let us show that p 1 > p 2 . 
(p 2 — p 1 ) satisfies 

j < D(p 2 - p 1 ), q > dtda + J o [e 2 (p 2 , q) - e^p 1 , q)]dtda = (21) 

(p 2 -p 1 )(0,a,/) = (22) 
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(p 2 -p 1 )(t,0,/) = (23) 

Let us choose q = (p 2 — p*) + . From the equality zfpf — z\p\ = z\{p\ — p\) + 
p1{zf - z\) follows that 

e1{v 2 M-p\) + )-e}{v\{pi-p}) + ) 

= c,((p 2 - p 1 ) + , (p1 -p\) + ) - q((p 2 - p 1 )-, (p? -p}) + ) 

+bi((pi - p\)\ (pi - p\) + ) - uiti - pir, (pi - p\) + ) 

+ ^piU!-flM-p]) + di 

We have already shown in the proof of Lemma 5.1 that Q((p 2 — p 1 ) - , (p 2 — 
Pi) + ) < and that bi((pf — p\)~ , (pf — p\) + ) = 0. Moreover since p 2 is non- 
negative the last term of equality (24) is nonnegative. Then we obtain that 

eV, (P 2 - pY) - eV, (p 2 - p 1 )*) > ^((p 2 - p 1 )*, (p 2 - P y ). 

Since (p 2 — p 1 ) satisfies (22) and (23) it also holds that 

/ <D(p 2 -p 1 ),(p 2 -p 1 )+>^a>0, 
Jo 

so that 

/ o e 1 ((p 2 -p 1 ) + ,(p 2 -p 1 n^a<0, 
and using the coercivity of e 1 we finally obtain p 2 < p 1 . 

Step 3: Let p > be given in V. We define two sequences (p 1,ri ) n >i and 
(p 2 '")n>i by (p 1 ' 1 = p, pW = JV ,n ) and (p 2 - 1 = p, p 2 ' n+1 = j^p 2 '™). 
From step 2 follows that p 1,2 > p 2,2 . 

In addition to p 1,3 = J^p 1,2 and p 2,3 = J^p 2,2 , let us define q 3 = J^p 1,2 . 
The inequality p 1,3 > q 3 follows from step 2, whereas p 2 ' 3 < q 3 follows from 
step 1. Therefore p 1,3 > p 2 ' 3 . An induction then shows that p 1 '" > p 2,ri , Vn > 
1 and since the sequences converge to the solution p 1 and p 2 of problem (P) 
associated with the vector of mortality rates z 1 and z 2 respectively, the proof 
is complete. □ 



6 Concluding remarks 



In this paper we have investigated a multi-region nonlinear age-size struc- 
tured fish population model. The model was formulated in a generic way so 
that it can be potentially used for various fish species. We formulated an ini- 
tial boundary- value problem and proved existence and uniqueness of a positive 
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weak solution. We also proved a comparison result which shows that the varia- 
tions in the mortality rate in each region have consequences on the population 
of fish in every regions. 

Other important problems need to be addressed now and are currently un- 
der progress. The first one concerns the numerical implementation of this 
model. In order to integrate numerically system (l)-(4) we use the character- 
istic method. Indeed this system can be viewed as a collection of systems of 
parabolic equations on the characteristic lines 

S = {{t + s, a + s); s G (0, s max (t , a ))}, 

where (to, a®) G {0} x (0, A) U (0,T) x {0}. Each of these systems is then 
integrated in time with an operator splitting method using the Lie formula 
([15], [16]). 

The second problem concerns the estimation of the different badly known 
parameters of the model (growth, mortality and migration rates) from the 
data available for fisheries and mentioned in the Introduction. In order to solve 
numerically this inverse problem, the implementation of a variational data 
assimilation method is under progress. The objective is to obtain a synthetic 
representation of the real system combining theoritical knowledge (the model) 
and experimental knowledge (the data). 
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